Due to their remote location within the Russian High Arctic, little is known about the mass balance of ice caps on Severnaya Zemlya now and in the past. Such information is critical, however, to building a global picture of the cryospheric response to climate change. This paper provides a numerical analysis of the climate and mass balance of the Vavilov Ice Cap on October Revolution Island. Mass balance model results are compared with available glaciological and climatological data. A reference climate was constructed at the location of Vavilov Station, representing average conditions for the periods 1974–1981 and 1985–1988. The site of the station has a mean annual temperature of −16.5°C, and an annual precipitation of 423 mm water equivalent. The mass balance model was calibrated to the measured mass balance, and tested against the time-dependent evolution of the englacial temperatures (to a depth of 15 m). The mass balance model was then converted to a distributed model for the entire Vavilov Ice Cap. Model results predict the spatial distribution of mass balance components over the ice cap. Processes involving refreezing of water are found to be critical to the ice cap's state of health. Superimposed ice makes up 40% of the total net accumulation, with the remaining 60% coming from firn that has been heavily densified by refreezing.
Introduction
Observations indicate a general retreat of the Earth's glaciers and ice caps over the last century attributed to a concurrent warming of the global climate. In order to quantify the current rate of glacier change and associated sea-level rise, and to forecast future changes, numerical models are required. Models have shown that the Earth's largest ice sheets react relatively slowly to climate change, making their study less appropriate to gauge the link between ice masses and climate than smaller glaciers and ice caps such as those which exist in the Russian Arctic. The aim of this paper is to quantify the components of mass balance on the Vavilov Ice Cap in Severnaya Zemlya. Such work is a necessary first step toward predicting the response of ice masses in this poorly understood region to future climate change.
THE ICE CAPS AND CLIMATE OF SEVERNAYA ZEMLYA
Severnaya Zemlya is the most easterly glacierized archipelago in the Russian High Arctic, consisting of four main islands with a total area of 36,800 km2, about half of which is presently ice-covered (Dowdeswell et al., in press) (Fig. 1). The archipelago is located between latitudes of 73° to 82°N and longitudes 90° to 110°E, and lies approximately 60 km north of the tip of the Taymyr Peninsula on the Russian mainland.
Due in part to its inaccessibility, relatively little is known about the climate of Severnaya Zemlya. Measurements in the ablation zone of the Vavilov Ice Cap, on October Revolution Island, record a mean July temperature of −0.7°C for the three-year period between 1974 and 1976, and an average of −2.7°C for the short summer from June to August (Bryazgin, 1981). Observations indicate that conditions are predominantly cloudy with frequent fog and that summer snowfall is not unusual. A continuous time series of temperature since 1933 and more fragmentary measurements of precipitation have been recorded at Fedorova Meteorological Station, located on the northern tip of the Taymyr Peninsular (77°43′N, 104°17′E). The mean annual temperature at Fedorova is about −15°C, with a July average temperature of 1.5°C (Dowdeswell et al., 1997). Mean annual precipitation measured at the meteorological station is about 0.20–0.25 m, compared with a mean value of 0.4 m water equivalent (w.e.) measured in the ablation zone of the Vavilov Ice Cap for the period 1980–1989 (Bol'shiyanov and Makeev, 1995). Mean annual temperatures recorded at Fedorova show great interannual variability, although the period from 1940 to the mid 1950s was markedly warmer than the following three decades (Dowdeswell et al., 1997).
The only published mass balance data for ice masses on Severnaya Zemlya are from the Vavilov Ice Cap. Snow survey studies provide balance data for 10 years during the periods 1975–1981 and 1986–1988 (Barkov et al., 1992). Measurements indicate a high annual variability ranging from −0.63 to 0.46 m w.e., with a mean and cumulative balance of −0.03 m a−1 and −0.29 m, respectively. Long-term observations also indicate that the southern margins of the Vavilov Ice Cap have advanced 150 to 450 m between 1952 and 1985 (Barkov et al., 1992).
Koryakin (1986) compared maps and aerial photographs of the Russian Arctic compiled in 1936–1938 and 1952 with satellite images taken in 1976 to examine changes in ice cover. Data from the period of observation were extrapolated to represent the whole of the 20th century. It is reported that glaciers on the three main archipelagos in the Russian Arctic have been shrinking since the early 1900s, although the rate of retreat declines with distance from the primary source of precipitation in the North Atlantic. Govorukha et al. (1987) claim a reduction of about 500 km2 in ice cover on the archipelago between 1931 and 1984 (on Bol'shevik Island, the ice margin has receded up to 2.5 km, and Kroshka Glacier on Pioneer Island and Morskoy Glacier on Komsomolets Island have completely disappeared).
THE VAVILOV ICE CAP AND VAVILOV STATION
The Vavilov Ice Cap has a relatively simple form with an area of 1771 km2 (Fig. 2; Dowdeswell et al., in press). The summit of the ice cap has an elevation of 708 m a.s.l. and is located in the southeastern half of a long and narrow plateau orientated in a southeast to northwest direction. All of the Vavilov Ice Cap terminates on land. The western margins descend to an elevation of about 50 m a.s.l. and are located less than 1 km from the Kara Sea. The northern part of the ice cap ends on higher ground with an elevation of about 250 m a.s.l., while the eastern and southern margins have an elevation between 100 and 200 m a.s.l. The total ice volume of the Vavilov Ice Cap is 567 km3 (Dowdeswell et al., in press). This is equivalent to about 1.4 mm of global sea level. Ice thickness generally increases from the margins to the geographic center of the ice cap where the maximum measured ice thickness reaches 595 m.
The Vavilov Meteorological Station (627 m a.s.l.) was established in 1974 and operated until 1988 as part of the Russian glaciological field base, located in the northern part of the Vavilov Ice Cap (Fig. 2). A suite of data is available for this location, including detailed measurements of climate, mass balance, near-surface ice density and stratigraphy, and englacial temperature.
STRUCTURE OF THE PAPER
The paper is divided into three parts. In the first part, a one-dimensional numerical glacier mass balance model is described. In the second part, a “reference climate” for Vavilov Station is constructed from meteorological data and used to describe the climatic regime of the Vavilov Ice Cap. In the third part, the model is converted to two dimensions and used in conjunction with the reference climate to quantify the mass balance of the entire Vavilov Ice Cap.
Mass Balance Model
The mass balance model has been adapted from the model of Greuell and Konzelmann (1994). Full details of the model are available in Bassford (2002), and a schematic representation is provided in Figure 3. The model comprises two parts, in which surface and englacial processes are computed.
SURFACE PROCESSES
The surface processes part of the model deals with the exchange of energy between the atmosphere and the glacier surface and the calculation of snowfall accumulation. Since measurements of the energy balance components are not available for the Vavilov Meteorological Station, energy fluxes have to be calculated using parameterizations developed on other glaciers. Two criteria were used to select the most suitable parameterizations from the various methods documented in the literature. First, the variables used in the parameterizations must be available from meteorological data measured at Vavilov Station, namely: daily values of air temperature, cloud amount and height, humidity, and precipitation. Second, the energy fluxes should be derived using parameterizations developed in situations that are analogous to the climatic regime of Severnaya Zemlya.
Short-wave radiation, divided into direct and diffuse components, is important in regions with a high mean cloudiness, such as Severnaya Zemlya. The partitioning between the two components depends linearly on cloud cover, so that for clear-sky conditions 85% of the solar radiation is direct, while for a total cloud cover this value is only 20%. Topographic shading is not considered in the model, given that the ice cap has no surrounding valley walls. The receipt of direct solar radiation is also influenced by the projection of the surface in the direction of the sun, which is modeled using equations provided by Walraven (1978).
Glacier energy balance studies have shown that short-wave radiation is the main contributor to melt energy (e.g., Braithwaite and Olesen, 1990; Arnold et al., 1996; Oerlemans, 2000). The surface albedo is an important parameter in energy balance models designed to calculate glacier melt as it controls the amount of absorbed short-wave radiation. Ground-based and satellite-derived measurements indicate large spatial and temporal variations in the surface albedo of valley glaciers and ice caps (Reijmer et al., 1999; Knap et al., 1999; Oerlemans et al., 1999). The spatial variation in albedo of ice caps in northern Severnaya Zemlya can be clearly seen on Landsat images acquired during late summer, in which a series of snow and ice facies is present, characterized by marked changes in albedo, ranging from bare ice at the margin, through a thin slush zone, to snow on the summit of the main ice caps (Dowdeswell et al., in press). The key factor influencing the evolution of surface albedo in this investigation is the complex metamorphosis of the snowpack during the melting season. Variations in the density of the surface should represent the changes in grain size and water content that affect albedo. Therefore, the model calculates surface albedo as a linear function of the density of the top 10 cm below the surface. Cloud amount is included in the parameterization to account for the slight increase in albedo with cloud cover, following Greuell and Konzelmann (1994).
Variations in incoming long-wave radiation are due mainly to changes in water vapor and temperature of the atmosphere and are calculated using the approach of Kimball et al. (1982). Outgoing long-wave radiation is calculated by making an assumption that snow and ice radiate as black bodies in the infrared part of the spectrum. The flux of radiation emitted from a black body depends only on its temperature. Therefore, the flux of outgoing long-wave radiation can be calculated from the glacier surface temperature.
Turbulent eddies transfer heat in the form of sensible heat (temperature) and latent heat (water vapor) between the atmosphere and the glacier surface. This transfer process is complex and depends on wind velocity, surface roughness, and the atmospheric stratification (Stull, 1997). A bulk transfer method is used to calculate the turbulent heat exchange (e.g., Oerlemans, 1992; van de Wal and Oerlemans, 1994; Fleming et al., 1997). This approach assumes the sensible and latent heat fluxes are proportional to the difference in temperature and humidity between a standard measuring height (2 m) and the glacier surface.
The final surface process to be accounted for is the accumulation of snow and ice. A common approach to determining the state of precipitation is to assume that precipitation falls in solid form when the air temperature is below a critical value, usually between 0 and 2°C. In this study, the threshold is set at 0°C; chosen so that the ratio between precipitation falling as rain or snow in the model is consistent with observations at the Vavilov Meteorological Station.
ENGLACIAL PROCESSES
Englacial temperature and density variations are associated with meltwater percolation and refreezing. Calculations are made on a vertical grid (in 1-D experiments the model has one vertical dimension representing a column of snow and ice, which is divided into a number of equally wide cells) extending from the snow/ice surface into the glacier to a depth of at least 20 m. This part of the model has been modified to calculate superimposed ice formation.
Changes in the temperature of cells are calculated using a thermodynamic equation (Bassford, 2002). The thermodynamic equation is valid for the upper 20 m of a glacier because the influence of advection and conduction parallel to the surface as well as energy produced by ice deformation are negligible (Greuell and Konzelmann, 1994). A zero heat flux is assumed at the lower boundary of the vertical grid.
Part of the short-wave radiation penetrates the surface into deeper layers of snow and ice. This may have a considerable effect on the surface energy balance and the englacial temperature. If such penetration is important, the surface will receive less energy for melting while the underlying layers are heated. The amount of absorption and scattering below the surface depends on the wavelength of the radiation and on the physical properties of snow and ice, in particular density. Greuell and Konzelmann (1994) model short-wave radiation penetration by assuming that all radiation with wavelengths greater than 0.8 μm (about 36% of the total short-wave radiation penetrating the glacier surface) is absorbed by the surface layer. The remaining part is extinguished at depth according to Beer's law (Bassford, 2002).
Surface melting and the percolation and refreezing of meltwater is treated in a manner similar to the model of Greuell and Konzelmann (1994), but includes an explicit calculation of superimposed ice formation. Melting is considered after calculating the energy exchange between the glacier surface and the atmosphere. If the temperature of the uppermost grid cell is raised above 0°C, then it is set to 0°C and the excess energy is used for surface melting. Subsurface melting is not calculated explicitly in the model but may occur as a result of radiation penetration. For conservation of energy, subsurface cells with a temperature >0°C are set to 0°C and the extra energy is used to heat the uppermost grid cell. If meltwater enters a cell with a temperature <0°C then refreezing takes place. The maximum refreezing rate that can occur in a cell at a point in time is limited by the amount of latent energy required to raise the cell temperature to 0°C. If the amount of meltwater percolating into a cell is less that its capacity, then all of the meltwater refreezes. If the quantity of meltwater exceeds the capacity then excess meltwater will be issued to the cell below. Refreezing results in an increase in cell density and temperature. If sufficient refreezing occurs in a cell to increase its density to that of ice (ρ = 910 kg m−3) then an impermeable ice layer forms in the snowpack. Meltwater continues to move down through the grid until either all of the meltwater refreezes or an impermeable ice cell is encountered, such as an ice layer within the snowpack or the underlying continuous glacier ice, on which superimposed ice forms (Wakahama et al., 1976). The impermeable ice surface is referred to as the snow-ice interface.
Greuell and Konzelmann (1994) use an empirical formula to calculate the densification of dry snow due to settling and packing. These processes are important in the absence of meltwater, such as in the dry snow zones of the Greenland and Antarctic ice sheets. Melting increases the densification rate of snow and firn. The refreezing of meltwater to form ice layers and superimposed ice represents the most rapid transformation of snow/firn to ice. For situations where meltwater refreezing is prevalent, it dominates the densification process (Paterson, 1994). This is generally the case in the accumulation area of ice caps on Severnaya Zemlya, which are characterized by a thin snow/firn pack, no more than 2–3 m deep, with frequent ice layers (Korotkevich et al., 1985). Therefore, densification of snow/firn occurs in the model by melting and refreezing while changes in density due to settling and packing are assumed to be negligible.
NUMERICAL DETAILS
Energy balance calculations are made on a staggered vertical grid so that fluxes are evaluated across the boundary between grid cells. Designing a grid system for the situation is non-trivial because of the moving upper boundary caused by snow accumulation and surface melting. The grid system consists of two integrated parts, one for the snowpack and a second for the underlying continuous glacier ice. The latter part extends to 20 m below the ice surface and consists of 25 cells of variable size ranging from 5 cm at the top to about 2 m at the lower boundary of the grid. The grid resolution is finest at the top of the grid because the gradients in ice temperature are greatest close to the ice surface. At greater depths below the ice surface, where the gradients in temperature are small, it is acceptable to have larger grid cell size.
The other part of the grid deals with the snowpack and has a regular cell size of 5 cm, except for the uppermost cell whose size is between 0 and 10 cm. As the depth of the snowpack rises and falls due to snow accumulation and melting, the number of grid cells above the ice surface changes. If a single snow cell remains, then the cell size is allowed to fall toward zero until all of the snow is melted.
Model simulations start at the beginning of the mass balance year, assumed to be “day of year” 274 (1 October). At the first time step of the model simulation, the snowpack is set to a depth of 5 cm with a density of fresh snow. At the end of the melting season, any snow remaining is removed to leave a bare ice surface in order to prevent the buildup of a deep firn layer. This approach is justified because under the present climate the actual depth of the snow/firn layer on the ice caps of interest is generally little more than the annual snowpack (Korotkevich et al., 1985).
A Reference Climate for Vavilov Station
The reference climate consists of an annual cycle of daily values for each of the climate variables required as input by the mass balance model, namely, air temperature, temperature range, precipitation, cloud cover, cloud base height, and relative humidity. Constructing a reference climate representative of average conditions at Vavilov Station for the periods 1974–1981 and 1985–1988 is complicated because the available data are fragmentary.
AVAILABLE CLIMATE DATA
Unpublished meteorological records for the Vavilov Station were retrieved from archives held at the Arctic and Antarctic Research Institute in St. Petersburg, Russia. The data set consists of daily values for various climatic parameters for the years 1978–1988. Longer running climate observations have been made from 1954 until the present at Golomyanniy Meteorological Station, located on the western tip of Sredniy Island (79°33′N, 90°38′E), 95 km northwest of Vavilov Station at an elevation of 7 m a.s.l, from which a continuous time series of mean monthly temperatures for the period 1920–1997 is available.
CONSTRUCTING THE REFERENCE CLIMATE
Temperature
Since meteorological data for Vavilov Station are available only for the years 1978–1988, temperatures were reconstructed for periods without data using the record at Golomyanniy as follows. Monthly mean air temperatures were calculated from daily values recorded at Vavilov Station. Following the approach of Lefauconnier and Hagen (1990), linear regression analysis was used to relate these monthly air temperatures to those measured at Golomyanniy Station. Despite a high value correlation coefficient (r = 0.99), the relationship between monthly temperature at the two climate stations varies depending on temperature. Therefore, in order to obtain a better fit between the two data sets, separate regression equations were calculated for each month of the year. Based on these regression equations, the mean monthly temperature was reconstructed at Vavilov Station for periods without measurements to form a complete time series of data from 1974–1981 and 1985–1988.
Precipitation
Since precipitation data recorded at Vavilov Station fail to represent interannual variations in winter mass balance over the ice cap, snowfall accumulation measurements were used to derive precipitation data for the reference climate. It is anticipated that measurements averaged over several years represent the seasonal distribution of precipitation.
A time series of monthly precipitation from September to May for the years 1974–1981 and 1985–1988 was obtained by multiplying the normalized average monthly precipitation by the winter mass balance. Precipitation for the summer months cannot be derived from the summer mass balance because it is difficult to separate the precipitation and melting components. Therefore, monthly precipitation from June to August is set to average values for Vavilov Station.
Other Variables
Monthly mean values of diurnal temperature range, cloud cover, cloud base height, and relative humidity were calculated from daily values measured at Vavilov Station. Because linear regression, as used to reconstruct temperature, failed to provide a reliable method of reconstructing data, periods without data were given a value equal to the mean of the monthly values.
Interpolating Daily Values from Monthly Means
Two final steps were taken to produce the reference climate. First, monthly values for the reference climate were obtained by calculating the average of all the monthly data. Second, a polynomial curve was used to interpolate daily values from the mean monthly data.
CLIMATE AT VAVILOV STATION
The reference climate at Vavilov Station can now be used to describe the climatic regime of the ice cap, together with daily measured data for 1987, which help to illustrate daily fluctuations concealed in the averaged data (Figs. 4 and 5). Conditions at the station are relatively cold and dry, with a mean annual air temperature of −16.5°C and an annual precipitation of 423 mm w.e. Temperatures fall to a minimum in February with a mean monthly value of −29.2°C, although daily temperatures as low as −40°C are common in winter (Figs. 4 and 5). A noticeable feature in Figure 5 is the large and abrupt fluctuations in winter temperature. The peaks in temperature generally correspond to precipitation events, with a concurrent increase in wind speed and cloudiness, as well as a fall in pressure. These conditions are indicative of cyclonic activity, originating in the North Atlantic where the large-scale circulation is dominated by a low pressure regime centered over Iceland. The frequency of these depressions is a key factor in determining the total winter snowfall in Severnaya Zemlya. Much colder and drier conditions prevail at other times on the Vavilov Ice Cap resulting from the high-pressure Siberian anticyclone that dominates the region during winter (Bryazgin, 1981).
In summer, air temperatures are generally highest in July when the mean monthly value rises to −0.5°C (Fig. 4). Between late June and August, daily temperatures fluctuate around the melting point, but cool periods often occur when temperatures fall well below 0°C (Fig. 5). Therefore, summer snowfalls are not uncommon. Conditions are generally cloudy with a mean cloud cover of about 0.82 and a high occurrence of fog (Fig. 4). The highest air temperatures and clearest skies occur when warm air masses heated over the mid-Siberian Plateau move northward over Severnaya Zemlya under anticyclonic conditions (Gordeichuk, 1997). Depressions from the North Atlantic continue to penetrate the region during summer causing cloudy and cool conditions, but lack the intensity of winter cyclones (Bryazgin, 1981).
Marked warming events lasting several days may take place in September, which interrupt the steady decrease in temperature after the end of summer (Fig. 5). During these warming events, surface melting may resume. Over 30% of the annual precipitation falls in September and October when the frequency of cyclones is highest (Figs. 4 and 5).
Calibration of the Mass Balance Model
Following the approach of Greuell and Konzelmann (1994), the performance of the 1-D mass balance model was evaluated by comparing model calculations with measurements of mass balance and ice temperature. The average observed mass balance at the site of the Vavilov Meteorological Station is 5 cm w.e. a−1 for the periods 1974–1981 and 1985–1988 (Barkov et al., 1992). At the same location, a time series of englacial temperature measurements was collected between December 1978 and September 1979 to a depth of 80 m below the ice cap surface (Barkov et al., 1988). Measurements indicate that ice temperature below a depth of 15 m varies within ± 0.1°C in response to seasonal changes in air temperature. The 15 m ice temperature, found to be −11.5°C, should thus be representative of at least a few years. Therefore, in addition to the average observed mass balance, the 15 m ice temperature is used as a diagnostic for testing the mass balance model.
An initial simulation was performed in which the 1-D mass balance model was run to steady state using the reference climate. For a first run, the model performed well, calculating a mass balance and 15 m ice temperature of −2 cm w.e. a−1 and −11.5°C, respectively. In order to improve the simulation of measurements, the model was calibrated by adjusting the albedo of fresh snow (details of the model's sensitivity to changes in other meteorological variables are provided in Bassford (2002)). The measured mass balance of 5 cm w.e. a−1 was obtained by increasing the albedo of fresh snow from 0.75 to 0.76, with the side effect of a slight decrease in 15 m ice temperature to −11.8°C (still a very close match to the measured value of −11.5°C).
In order to evaluate the ability of the calibrated model to simulate the evolution of englacial temperature, the model was tested against a time series of ice temperature measurements obtained close to Vavilov Station between December 1978 and September 1979. During this period, measurements were taken from a string of thermistors installed at 22 levels with a variable spacing of 1, 5, and 10 m at depths of 0–15, 15–30, and 30–80 m, respectively (Barkov et al., 1988). For comparison, the evolution of englacial temperature was simulated by the mass balance model using daily climate data for the mass balance year 1978–1979. Calculated ice temperature data, corresponding to the same dates as measured profiles, were output from the model (Fig. 6). Although there are some noticeable discrepancies, particularly in late summer, the calculated ice temperatures show the same general evolution as the measurements. Both data sets have the same ice temperature range at the surface and a similar decrease in the magnitude of variation with depth (Fig. 6).
The Annual Cycle at Vavilov Station
The calibrated 1-D mass balance model can be used with the reference climate to simulate an annual cycle of the surface mass balance at Vavilov Station. The averaging involved in constructing the reference climate means that the extremes of certain variables, such as daily melt rates, will be underestimated. However, this problem should not invalidate results that are integrated over longer time periods, e.g., weeks to months. Although in reality these cycles vary from year to year, this simulation helps to identify and explain the general features of the annual cycle.
Melting starts in early June and ends in late August, on “days of year” 159 and 241, respectively (Fig. 7). The peak melting rate of only 0.91 cm w.e. day−1 takes place on “day of year” 200, coinciding with the maximum mean daily temperature of 0°C. In reality, the maximum melt rates are likely to be much higher than this on days with a mean temperature several degrees above zero. Note that melting occurs even when the mean daily temperature is below 0°C. This is because the surface energy balance is positive and thus provides energy for melting despite subzero temperatures. Diurnal variation also means that air temperatures may rise above 0°C even though the mean daily value is below the melting point. The rapid decrease in melt rates after the mid-summer maximum is interrupted by a second peak (Fig. 7). This feature occurs because of a decrease in surface albedo which temporarily increases melt rates.
Runoff is delayed by 21 days after the onset of melting because of refreezing in the snowpack (Fig. 7c). Even then, the amount of runoff is considerably less than melting because of continued refreezing in the snowpack and superimposed ice formation. Night time cooling means that refreezing in the snowpack continues for much of the summer, although the amount of refreezing decreases rapidly after the winter cold content of the snowpack is eliminated. Melting and refreezing in the snowpack is reflected by a steady decrease in snowpack depth with a synchronous increase in mean density (Fig. 7d). Initially, all meltwater percolating to the ice surface refreezes as superimposed ice, but the rate of refreezing on the ice surface quickly decreases to a negligible amount as the ice surface temperature is raised to 0°C by the release of latent heat (Fig. 7). For a brief period between “days of year” 210 and 220, all meltwater produces runoff, the amount of which is augmented by precipitation falling as rain. The snowpack melts away completely by “day of year” 230, but 10 days later snow starts to accumulate again. A period of renewed superimposed ice formation lasting about two weeks occurs toward the end of the melt season as the ice surface temperature drops below 0°C, although the rate of refreezing is very small.
The total amount of surface melting during the balance year is 51.3 cm w.e. a−1, with rainfall contributing an additional 1.2 cm w.e. a−1 of water to the system. About 70% or 37.2 cm w.e. a−1 of the combined total of meltwater and rainfall produces runoff. The difference is due to 8.5 cm w.e. a−1 of meltwater refreezing in the snowpack and 6.8 cm w.e. a−1 of superimposed ice formation. This is consistent with stratigraphic data from a snow pit and a shallow ice core at the site of the Vavilov Station, indicating layers of heavily densified firn and superimposed ice with an average thickness of about 10 cm (Korotkevich et al., 1985).
The calculated amounts of refreezing in the snowpack and superimposed ice formation may seem small but are significant given the relatively low mass turnover at Vavilov Station. In fact, the net accumulation of 5 cm w.e. a−1 is made up entirely of superimposed ice.
The Distributed Mass Balance Model
The 1-D mass balance model was transformed into a 3-D distributed model following Arnold et al. (1996), allowing the mass balance of the whole ice cap to be calculated.
LAPSE RATE
A temperature lapse rate is used in the mass balance model to determine the value of air temperature with altitude. It is an important parameter in the model because several processes are temperature dependent, such as the fluxes of sensible heat and incoming long-wave radiation, as well as the form of precipitation. Ideally, the temperature lapse rate of the free atmosphere should be measured by repeated radiosonde balloon soundings, but reliable data of this kind were unavailable for the study area. Therefore, the temperature lapse was set to a standard value of 0.006°C m−1.
CALIBRATION
Russian measurements of winter snow accumulation along two mass balance profiles indicate that the southern part of the ice cap receives about 50% more precipitation than north-facing slopes. However, there are no measured data available to establish the gradients in precipitation across the eastern and western parts of the ice cap. Therefore, the spatial variation in precipitation was determined by treating precipitation as a tuning variable in the model which was varied to improve the match between the calculated and measured mass balance. Measurements presented in Barkov et al. (1992) were used for this purpose, consisting of a map of the distribution of mean net mass balance over the entire Vavilov Ice Cap for the period 1974–1988 (Fig. 8).
Calibration involved running the distributed mass balance model with the reference climate for 440 points spaced 2 km apart over the ice cap. At each point, precipitation was initially set to an amount equal to the reference climate. Then, if necessary, the model was tuned to calculate the measured specific mass balance to within ±1 cm w.e. a−1 by adjusting the amount of precipitation. The temporal distribution of precipitation over the ice cap was assumed to be the same as at the Vavilov Station, although the form of precipitation was varied according to air temperature. This procedure was performed automatically using an iterative process called the bisection method (Press et al., 1989). A bicubic spline was then used to interpolate between the 440 points to produce a value of annual precipitation for each cell of a 2-D grid over the ice cap.
According to the model, there is a significant gradient in precipitation across the Vavilov Ice Cap, with annual precipitation decreasing from a maximum of 63 cm w.e. in the southwest to a minimum of 31 cm w.e. in the north compared to an average value of 49 cm w.e. (Fig. 9). This gradient in precipitation is evident from the altitude of the transient snowline identified in a late summer Landsat image, which varies from about 350 m a.s.l. in the southern part of the ice cap to over 650 m a.s.l. in the north (Bassford, 2002). In many mass balance modeling studies the variation in precipitation over the surface of a glacier is expressed as a simple function of altitude (e.g., Oerlemans, 1992; Fleming et al., 1997), but in the case of the Vavilov Ice Cap altitude appears to have a second-order effect on the distribution of precipitation. The greater amount of precipitation in the southwestern part of the ice cap is explained by the direction of prevailing winds, which bring moist air masses from the southwest, associated with depressions originating in the North Atlantic. The Vavilov Ice Cap is the first topographic obstacle that moist air masses encounter as they pass over the archipelago. Therefore, the southern half of the ice cap receives relatively large snowfalls while the north and eastern parts receive considerably less snowfall because they are located in a precipitation shadow. Local anomalies to this general distribution of precipitation occur as a result of the redistribution of snow by wind. For example, the large amount of annual precipitation close to the southwestern margin is likely to occur because blown snow accumulates in a hollow formed between the south- and southwest-facing slopes of the ice cap (Fig. 9).
THE REFERENCE STATE
The calibrated distributed model was run with the reference climate to calculate the net mass balance of the Vavilov Ice Cap. The modeled mean net balance of the entire ice cap is −2.2 cm w.e. a−1, which compares closely with a measured average value of −2.8 cm w.e. a−1, indicating that the ice cap was approximately in balance for the period 1974–1988 (Barkov et al., 1992). Net mass balance calculated by the model varies from −58 cm w.e. a−1 at 42 m a.s.l. on the western margin of the ice cap to 35 cm w.e. a−1 at a point located approximately 3 km southwest of the summit with an altitude of 685 m a.s.l. (Fig. 10). The mean equilibrium line altitude is 498 m a.s.l. but ranges from 350 to 500 m a.s.l. in the southern half of the ice cap to a maximum of 621 m a.s.l. in the north, reflecting the marked gradient in precipitation across the ice cap (Fig. 10a). In agreement with Barkov et al. (1992), the accumulation zone covers a planimetric area of about 810 km2 which is 45% of the total area of the ice cap.
THE IMPORTANCE OF MELTWATER REFREEZING
Components of the mean net balance calculated by the model are listed in Table 1 and distributed results are shown in Figure 10. Surface melting varies from just over 100 cm w.e. a−1 at the lowest elevation of the ice cap on the western margin to about 45 cm w.e. a−1 close to the summit. However, the amount of surface melting is not a simple function of elevation. A greater amount of surface melting occurs in the northern part of the ice cap relative to the same elevation in the south because melting is reduced in the south by the higher surface albedo associated with a longer lasting summer snowpack.
On average, 81% of meltwater is lost from the ice cap as runoff, with the remainder refreezing in the snowpack or on the ice surface as superimposed ice (Table 1, Figs. 10c and 10d). Both of these processes are significant over the whole ice cap surface (Figs. 10e and 10f). The amount of refreezing in the snowpack increases from 5 cm w.e. a−1 at the margins of the ice cap to a maximum of 11 cm w.e. a−1 at the summit, because the cold content of the snowpack at the start of summer is greater at higher elevations where the snowpack is deeper and the air temperature is lower. Figures 10a and 10f indicate that the amount of superimposed ice formation, which varies from 6 to 10 cm w.e. a−1, is closely associated with winter snow accumulation. As discussed earlier, the formation of superimposed ice in the model can only take place while a snowpack is present. Therefore, the maximum amount of superimposed ice formation occurs in the southeastern part of the ice cap where snow accumulation is greatest.
Although all of the meltwater and rainfall that refreezes in the ablation area is subsequently melted and lost from the ice cap as runoff, refreezing is still an important process in this part of the ice cap since it delays the start of runoff and so reduces net ablation. Above the equilibrium line, superimposed ice formation contributes 40% of the total net accumulation of the ice cap, with the remaining 60% comprising of firn that has been heavily densified by refreezing (Table 1). Model results indicate that the superimposed ice zone covers an area of 179 km2, approximately 20% of the accumulation area (Fig. 11). The extent of the superimposed ice zone is reflected by a difference of 36 m between the mean altitude of the equilibrium and snow/firn lines, located at 498 and 534 m a.s.l., respectively.
The large contribution of refreezing to net accumulation on the Vavilov Ice Cap is also evident from stratigraphic analysis of an ice core extracted from the accumulation zone (Fig. 11). Korotkevich et al. (1985) identified three different types of ice in the core based on grain size and the proportion of air inclusions in the ice which reflect the process by which the ice formed. In the upper 5 m of the core, 24% of the section was found to be regelation ice (formed by the compression of dry firn), 43% was identified as infiltration ice (formed by the meltwater refreezing at the start of the melting season), and 33% of the core was infiltration-congelation ice (formed by the refreezing of saturated firn in late summer).
An important implication of these results is that care must be taken when estimating the mass balance of Arctic glaciers and ice caps using satellite images such as those acquired by Landsat. Misinterpretation of the snow/firn line as the equilibrium line is likely to result in an underestimation of the mass balance of an ice mass. Clearly, superimposed ice formation and meltwater refreezing in the snowpack are key processes which should be incorporated by mass balance models applied to Arctic ice masses. Neglect of these processes in the model used in this study would result in a decrease in the calculated mean net balance from −2.2 cm w.e. a−1 to −18.0 cm w.e. a−1, a significant change given the small mass turnover of the Vavilov Ice Cap.
Summary and Conclusions
A reference climate was constructed for Vavilov Meteorological Station, representing average conditions for the periods 1974–1981 and 1985–1988. The Vavilov Station has a mean annual temperature of −16.5°C and an annual precipitation of 423 mm w.e. Winter conditions are cold and dry with infrequent snowfall events associated with cyclonic activity, while summers are cool and cloudy with temperatures fluctuating around 0°C. Temperature inversions caused by a surface radiation deficit are common during winter, particularly in March when temperatures at Vavilov Station (627 m a.s.l.) are on average 2.3°C higher than at Golomyanniy Station (7 m a.s.l.).
A 1-D mass balance model was successfully calibrated to calculate the annual average mass balance (5 cm w.e. a−1) and 15 m ice temperature (∼11.5°C) at Vavilov Station. The calibrated model proved to be well capable of simulating the measured evolution of englacial temperature in the upper 15 m of ice between December 1978 and September 1979.
The calibrated 1-D model was used to simulate an annual cycle of the surface energy balance, englacial temperature, and mass balance for the reference climate at Vavilov Station. The energy balance in winter is characterized by a radiation deficit balanced by a sensible and subsurface heat flux toward the surface. Net radiation is positive from April to September and is balanced by heating and melting of the surface, together with a sensible heat flux directed toward the atmosphere. In summer, the decrease in albedo associated with surface melting and densification of the snowpack is the most important factor determining the radiation balance. Seasonal variation in the surface energy balance affects englacial temperature to a depth of about 15 m where the calculated temperature remains constant at −11.8°C, about 5°C warmer than the mean annual air temperature. Model results indicate that the key factor responsible for this difference is the insulation provided by the winter snowpack, while a large amount of heating in summer is caused by meltwater refreezing and short-wave radiation penetration and absorption.
The mass balance model calculates that about 70% of meltwater produces runoff, with the remainder refreezing in the snowpack and on the ice surface as superimposed ice. Therefore, meltwater refreezing is an important component of the mass balance at Vavilov Station, and in fact the annual net accumulation of 5 cm w.e. is made up entirely of superimposed ice. Sensitivity experiments show that snow depth at the onset of melting and summer air temperature are key factors controlling the quantity of meltwater refreezing. However, a significant amount of refreezing occurs in a wide range of conditions and is therefore likely to be important in most parts of the Vavilov Ice Cap in the majority of balance years.
A 3-D distributed mass balance model was developed from the 1-D model and calibrated for the Vavilov Ice Cap by treating precipitation as a tuning variable. The calibrated model was used to simulate the mass balance of the Vavilov Ice Cap for the reference climate defined at Vavilov Station. In close agreement with field measurements, the model calculates that the mean net mass balance of the ice cap is −2 cm w.e. a−1 for the period 1974–1988, with a mean equilibrium line altitude of 498 m a.s.l. The model confirms that meltwater refreezing in the snowpack and superimposed ice formation are important processes over the entire ice cap, reducing the amount of meltwater producing runoff by an average of 19%. Superimposed ice is estimated to contribute 40% of the total net accumulation of mass on the ice cap, with the remaining 60% comprising comprised of firn that has been heavily densified by refreezing.
The mean net balance of the Vavilov Ice Cap shows a high interannual variability which is caused primarily by variations in the amount of summer melting. Model results indicate a large year to year variation in the equilibrium line altitude, which sometimes exceeds the summit of the ice cap. Intense surface melting in the accumulation zone during warm summers prevents the buildup of a thick firn layer by rapid transformation of firn to ice through refreezing and by removing mass through runoff.
Acknowledgments
We thank Lisette Klok, an anonymous referee, and the editor for providing very helpful reviews. Bassford acknowledges funding from a U.K. NERC Ph.D. studentship. Dowdeswell acknowledges funding from NERC grant GR3/9958 and the EU FP6 SPICE Project.
References Cited
FIGURE 1.
Map of Severnaya Zemlya showing the main ice caps in the archipelago. Inset is the location of Severnaya Zemlya within the Eurasian High Arctic and the nearby Russian Arctic archipelagos of Franz Josef Land and Novaya Zemlya
FIGURE 2.
The surface topography of the Vavilov Ice Cap. Contours are at 50 m intervals. The location of Vavilov Station is denoted by a “×” in (a). Adapted from Dowdeswell et al. (2002)
FIGURE 3.
Schematic structure of the surface mass balance model. M denotes transient mass balance, P* is the precipitation rate, and R is runoff
FIGURE 6.
Measured and modeled ice temperature profiles at the Vavilov Station for the period December 1978 to September 1979. Measured data are from Barkov et al. (1988)
FIGURE 7.
Model results for the reference climate at the Vavilov Station. (a) Daily mean temperature of the air (screen level), snow/ice surface and ice surface. (b) Depth and mean density of the snowpack. (c) Surface melting and runoff. (d) Meltwater refreezing in the snowpack and superimposed ice formation
FIGURE 8.
Specific mean mass balance (cm w.e. a−1) of the Vavilov Ice Cap for the period 1974–1988. The black triangle marks the location of Vavilov Station, the accumulation area is shaded by diagonal lines, and bare land is stippled. The figure was redrawn from a map constructed by Barkov et al. (1992) based on mass balance measurements along four profiles, labeled 1–8, and aerial photographs. Measurements along each profile were conducted for different time periods detailed in the table below. The only available data from these measurements are the mean winter and net mass balance at points along profiles 1–2 and 3–4 for the period 1974–1981
FIGURE 9.
Annual precipitation over the Vavilov Ice Cap, determined by tuning the mass balance model. Isolines are at 5 cm w.e. intervals. Bare land is stippled
FIGURE 10.
(a) Modeled mass balance of the Vavilov Ice Cap under the reference climate. Contours are at 10 cm w.e. a−1 intervals. (b) Modeled snowfall on the Vavilov Ice Cap under the reference climate. Contours are at 5 cm w.e. a−1 intervals. (c) Modeled surface melting of the Vavilov Ice Cap under the reference climate. Contours are at 10 cm w.e. a−1 intervals. (d). Modeled meltwater runoff from the Vavilov Ice Cap under the reference climate. Contours are at 10 cm w.e. a−1 intervals. (e) Modeled refreezing in the snowpack on the Vavilov Ice Cap under the reference climate. Contours are at 1 cm w.e. a−1 intervals. (f) Modeled superimposed ice formation on the Vavilov Ice Cap under the reference climate. Contours are at 1 cm w.e. a−1 intervals
